Two-Dimensional Dynamics of Ultracold Atoms in Optical Lattices 
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We analyze the dynamics of ultracold atoms in optical lattices induced by a sudden shift of the 
underlying harmonic trapping potential. In order to study the effect of strong interactions, dimen- 
sionality and lattice topology on transport properties, we consider bosonic atoms with arbitrarily 
strong repulsive interactions, on a two-dimensional square lattice and a hexagonal lattice. On the 
square lattice we find insulating behavior for weakly interacting atoms and slow relaxation for strong 
interactions, even when a Mott plateau is present, which in one dimension blocks the dynamics. On 
the hexagonal lattice the center of mass relaxes to the new equilibrium for any interaction strength. 
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The possibility to confine ultracold atoms in optical 
lattices has opened up a new research area, where inter- 
acting quantum many-body systems consisting of bosons 
[H, 0] and fermions [3, 0] can be studied with unprece- 
dented high precision and tunability Q. A new and ex- 
citing development is to study the dynamics and out-of- 
equilibrium behavior of those systems, see e.g. 0, which 
can also be used as an experimental probing technique 
0, 0|- In particular, it is possible to bring the system far 
out of equilibrium by performing an instantaneous shift 
in the underlying harmonic trapping potential. In this 
way particle transport can be investigated. This proce- 
dure is schematically illustrated in Fig. [T] Experiments 
in this setup with a three-dimensional Bose- Einstein con- 
densate subject to a one-dimensional optical lattice re- 
vealed a transition from coherent to dissipative behavior 
El, [i2|- For a truly one-dimensional gas, the im- 
pact of the optical lattice was found to be even stronger: 
already for small lattice depths dissipative motion was 
observed 13, which has been theoretic ally investi- 
gated with a variety of methods [H, [B, [13, [ll, [H, [l^l . 
A three-dimensional fermionic cloud subject to a one- 
dimensional optical lattice showed strongly suppressed 
center of mass motion 
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24l | and experimental 




22| . Related theoretical 
25| work on a Bose gas in 
a moving optical lattice showed a momentum-dependent 
breakdown of superfluid motion. 

However, both experimental and theoretical research 
in this direction has up to now [iit been restricted to a 
one-dimensional optical lattice. In this Letter we ana- 
lyze the behavior of repulsively interacting bosons in two 
dimensions, where geometrical considerations play a pro- 
found role. The influence of the lattice dimensionality is 
most clearly seen when a Mott plateau is formed. Unlike 
in one spatial dimension, where this leads to insulating 
behavior [l6l |. we find that in two dimensions there is 
always relaxation of the center of mass to the new equi- 
librium position. This behavior even persists in the limit 
of hard-core bosons. Moreover, in two dimensions, one 
can choose different lattice geometries, which severely in- 
fluences the dynamics. Experimentally this is possible by 
varying the angle between the laser beams which make up 



FIG. f : (Color online) Schematic of the setup: the harmonic 
potential is suddenly shifted, bringing the system far out of 
equilibrium. In the insets the two lattice geometries we study 
are shown: the square lattice (left) and the hexagonal lattice 
(right). 



the optical lattice. Here we compare the behavior on the 
square lattice to the hexagonal lattice and find remark- 
ably different behavior. This originates from the single 
particle dynamics, which on the square lattice can be de- 
scribed in terms of Bloch oscillations. Those are absent 
on the hexagonal lattice, where the atoms move along 
the cquipotcntial lines. Therefore, strongly interacting 
bosons on the hexagonal lattice show a much quicker re- 
laxation of the center of mass than on the square lattice. 
In order to study bosonic atoms with arbitrarily strong 
repulsive interactions we apply the Gutzwiller mean-field 
theory. This includes the limit of infinite interaction 
strength (hard-core bosons). 

For a deep optical lattice and moderate filling, the 
bosons can be described by the single band Bose- 
Hubbard Hamiltonian 
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where bi is the bosonic annihilation operator at site i, 
/i is the chemical potential, J is the hopping amplitude 
and U is the on-site repulsion. J and U can be expressed 
in terms of the atomic interparticle scattering length a, 
mass m and laser wavelength and intensity We will 
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use them as effective parameters. V{i,t) is the underlying 
harmonic potential which we take equal to: 



y(z,t) = l/o|x.-xo(i)|2 



(1) 



In the following we set J = 1 and use U = U/J,fi~i-i/J 
and Vo = Vq/J as dimensionlcss parameters. Time is 
expressed in units of h/ J. The shift A in the harmonic 
potential is expressed in terms of the lattice spacing a. 

For weakly interacting bosons, a Bose-Einstein con- 
densate forms, whose dynamics can be described by the 
Gross-Pitaevskii equation. A complementary approach, 
which is also valid for strongly interacting bosons, is the 
time-dependent Gutzwiller technique [27| where the cou- 
pling between the lattice sites is treated in a mean-field 
approximation. For inhomogeneous systems, this proce- 
dure has to be carried out in a space-resolved version, 
where with each site a separate order parameter is as- 
sociated. The total many-body wavefunction is within 

this approximation given as j'l') = Hi X^ij 



In practice, the infinite sum over the particle numbers 
n is replaced by a finite sum, by introducting a cut-off 
Nc depending on the strength of the interaction and the 
local density. The dynamics is g overned by the set of 
coupled differential equations 27 1 



m 

+ (^^n{n^l) + V{z,t)-f,^rn, (2) 

where = (6,) = En V^Un-i)* fl This Gutzwiller 
approximation is a highly efficient method for studying 
dynamics in higher dimensional lattices. It conserves en- 
ergy exactly and particle number with a very good accu- 
racy. The latter, however, is only true if the sites are se- 
quentially updated; a parallel update of all sites together 
violates particle number conservation. The validity of 
the Gutzwiller approximation is further justified by the 
fact that for small interactions it incorporates the Gross- 
Pitaevskii dynamics . Also the limit of strong inter- 
actions in one dimension is correctly reproduced by the 
Gutzwiller approximation (see Fig. [3|): in that case the 
dynamics of the center of mass is completely blocked [3l ■ 
The Gutzwiller approximation is naturally restricted to 
zero temperature, since it neglects phase fluctuations. 
We therefore only consider T = here. 

We first investigate the influence of the lattice topology 
on the single particle dynamics, which will reflect itself 
in the behavior of weakly interacting bosons. For small 
shifts, the bosons perform dipole oscillations around the 
shifted center. Interactions between the bosons lead to 
damping that increases with the strength of the interac- 
tion. This behavior is limited to small shifts, for which 
the shifted wavefunction still has an overlap with the 
groundstate-wavefunction in the shifted potential. Al- 
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FIG. 2: (Color online) Center of mass (COM) dynamics for 
non-interacting bosons after an instantaneous shift in the har- 
monic potential on a two-dimensional square lattice (left col- 
umn) and a hexagonal lattice (right column) obtained with 
exact diagonalization. For all plots Vo — 0.3. Shown are shifts 
A ^ 2 (red), A = 4 (green), A = 6 (blue), A = 8 (purple), 
A = 10 (light-blue) and A = 12 (yellow). In b) the long-term 
dynamics is depicted in the inset. In the two lower panels 
the density distribution of non-interacting bosons is plotted 
c) for the square lattice, (t/ J — 300) and d) for the hexagonal 
lattice {t/ J = 3000) for a shift of ten lattice sites. Note that 
for the hexagonal lattice the color coding is on a logarithmic 
scale. 



ternatively, one can argue that for larger shifts the po- 
tential energy introduced into the system by performing 
the shift cannot totally be converted into kinetic energy, 
which is restricted because of the single-band descrip- 
tion. For larger shifts the lattice structure becomes very 
important, which is most clearly seen in the limit of non- 
interacting bosons (see Fig. [2]). If J7 = 0, the Hamilto- 
nian on the square lattice is the sum of two commuting 
one-dimensional Hamiltonians: 'Hu=o ~ 'Hx + Ti-y with 
\Hx:'Hy\ = 0. This implies that the dynamics is effec- 
tively one-dimensional and the single particle eigenstates 
are products of the one-dimensional eigenstates, which 
are highly localized [2^. As a result, the atoms perform 
Bloch oscillations around the shifted position instead of 
dipole oscillations after a large shift [30[ as shown in Fig. 
[5^). Because of the non-linearity of the harmonic poten- 
tial, the Bloch oscillations contain multiple frequencies. 
Because the eigenstates are localized, the wavepacket re- 
mains localized as well (Fig. [5J:). The Bloch oscillations 
persist for small interactions between the bosons, which 
however cause the Bloch oscillations to be damped, due 
to dephasing [2^. For stronger interactions the Bloch 
oscillations disappear. Instead, the center of mass shows 
dissipative dynamics. This is possible, because in this 
case the potential energy can be converted into interac- 
tion energy, and the center of mass can relax to zero. 

The description in terms of Bloch oscillations does 
not apply to other lattice structures, where the Hamil- 
tonian cannot be decomposed into two commuting one- 
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FIG. 3: (Color online) Density distribution (color coding) 
and center of mass (green line) as a function of time for a 
one- dimensional system. Parameters are chosen as Vb = 0.3, 
A = 2, N = 20 and from left to right U = 0.2, [/ = 2 and 
U = 20. For the latter case a Mott plateau forms, which 
completely blocks the relaxation of the center of mass. 

dimensional parts. In particular, we have investigated 
the behavior of free and weakly interacting bosons on 
a hexagonal lattice. In this case we find for arbitrarily 
large shifts a relaxation of the center of mass to the min- 
imum of the harmonic potential (Fig. [Sb)). Inspection 
of the density profile shows that this occurs, because on 
the hexagonal lattice the atoms move along the equipo- 
tential lines. The reason for this behavior is, that on 
the hexagonal lattice plus harmonic potential the single- 
particle eigenstates form a ring-like structure. Therefore, 
the bosons perform a ring-like expansion after a large 
shift and although the center of mass relaxes to the new 
equilibrium, the atoms actually never do. 

It is worth noticing that spinless fermions behave in the 
same way as weakly interacting bosons. On the square 
lattice they perform Bloch oscillations after a large shift 
[2^ , whereas on the hexagonal lattice the center of mass 
relaxes to zero, whereas the density forms a ring-like 
structure. This means that, especially on the square lat- 
tice, spinless fermions behave qualitatively different from 
hard-core bosons. The latter show relaxation of the cen- 
ter of mass to the new equilibrium, as we will demon- 
strate later. This is unlike the one-dimensional situation, 
where hard-core bosons behave as free fermions. 

We now turn to strongly interacting bosons, for which 
the dimensionality of the lattice is very important. In one 
spatial dimension, the motion gets completely blocked 
when a Mott plateau is present, because the atoms can 
not pass each other and the center of mass never reaches 
the new equilibrium position [l6t . This behavior is re- 
produced by our Gutzwiller calculations as shown in Fig. 
[31 In two dimensions, however, the center of mass mo- 
tion shows dissipative dynamics and fully relaxes to the 
new equilibrium (Fig. 3]). We did extensive calculations 
for a wide range of parameters and always found this be- 
havior. The reason is that the superfluid shell can freely 
move around the Mott-insulating core, thus allowing the 
system to relax. This is visible in the density plots in 
Fig. [5l which show that initially the Mott plateau re- 
mains inert, but afterwards shows dynamical melting in- 
duced by the dynamics of the superfluid ring. Moreover, 
the insulating plateau scatters the atoms away from their 
single-particle trajectories, which leads to a significant 
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FIG. 4: (Color online) Center of mass dynamics as a func- 
tion of time for strongly interacting bosons on the square lat- 
tice (left column) and hexagonal lattice (right column) for 
iV = 100 and Vo = 0.3. The upper plots are for U = 40, 
where a Mott plateau is present before the shift. Displayed 
is the dynamics after a shift of ^ = 2 (red), A = 4 (green), 
A = 6 (blue), A = 8 (purple) and A = 10 (light-blue). The 
lower panels show the long-time dynamics for a large shift 
of A = 10 for different interaction strengths U — 2 (red), 
U =10 (green), U = 20 (blue), U = 40 (purple) and [/ = cx) 
(light-blue). In the inset of c) the long-time dynamics for the 
center of mass of hard-core bosons is shown on a logarithmic 
time-scale. In d) the U = 0-result is depicted in yellow for 
comparison. 
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FIG. 5: (Color online) Density distributions for a system of 
strongly interacting bosons after a shift in the harmonic po- 
tential on a square lattice with Vo = 0.016, A = 10 and 
A*' = 2551 and U = AO For these parameters, a very large 
Mott plateau forms in the center. Still the center of mass 
fully relaxes to zero as shown in e). The four panels display 
the density profile for a) t/J = 15, b) t/J = 30, c) t/J = 45 
and d) t/J = 9000. 



broadening of the density profile. Therefore, the final 
state is not the equilibrium state in the shifted poten- 
tial, because the total energy is conserved. In particular, 
the Mott plateau has disappeared (Fig. HJi)). The sup- 
pressed density in the center and the broadening of the 
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FIG. 6: (Color online) Radial density profile before the shift 
of the harmonic potential (red) and after relaxation to the 
new equilibrium (green) for the square lattice (left pictures) 
and hexagonal lattice (right pictures) for (7 = 2 (upper row) 
and U — oo (lower row). The blue line is a thermal density 
profile. The bosons on the square lattice have relaxed from a 
shift of six lattice sites; the bosons on the hexagonal lattice 
relaxed from a shift of ten lattice sites. Other parameters are 
chosen as iV = 100 and V = 0.3. 



tions, dimensionality and lattice topology on the particle 
transport of interacting bosons induced by a change in 
the underlying harmonic potential. We find that on the 
square lattice the single-particle dynamics on the square 
lattice is insulating after a large shift, and can be de- 
scribed in terms of Bloch oscillations. On the hexago- 
nal lattice the single particle dynamics shows relaxation 
along the equipotential lines. For strongly interacting 
bosons, we find that in two dimensions the presence of a 
Mott plateau does not prohibit relaxation of the center 
of mass to the new equilibrium, in contrast to the one- 
dimensional situation. This relaxation is very slow on the 
square lattice, but much faster on the hexagonal lattice. 
Our numerical predictions can be directly verified exper- 
imentally by performing absorption measurements, from 
which both the center of mass dynamics and the density 
profiles can be obtained. 

We thank Michael Kohl and Henning Moritz for useful 
discussions. This work was supported by the SFB-TRR 
49 of the German Science Foundation (DPG). 



density profile is shown in Fig. [B] We compare this with 
a finite temperature density profile, which is calculated 
neglecting phase fluctuations. We derive the temperature 
from the energy induced in the system by the shift in the 
potential. This leads to a good agreement, which shows 
that the density profiles for those parameters thermalize. 
The complete relaxation of the center of mass even per- 
sists for hard-core bosons. It is worth noticing that the 
time-scale for the relaxation cannot easily be expressed in 
terms of the interaction strength, because the relaxation 
is not due to single-particle tunneling processes. There- 
fore, the observed relaxation of the center of mass is very 
slow. This is partly due to the parameters chosen for 
the simulations; taking a more shallow harmonic poten- 
tial and a larger number of particles, as experimentally 
more realistic, brings the timescale back into the exper- 
imentally relevant regime. On the hexagonal lattice, the 
qualitative behavior of strongly interacting bosons is the 
same as above: the center of mass always relaxes to zero. 
However, quantitatively the lattice topology is very im- 
portant: the dynamics on the hexagonal lattice is much 
faster (Fig. IHd) and d)) and equal to the single-particle 
relaxation time. This is because the single particle dy- 
namics is along the equipotential lines. Therefore the 
atoms avoid each other and are able to relax faster. From 
Fig. [n]b) and d) it is obvious that this leads to a density- 
profile that is suppressed in the center and has the same 
ring-like structure as the non-interacting bosons obtain 
after relaxation of the center of mass. As an important 
consequence, this implies that the density profile does not 
thermalize on the hexagonal lattice, because the thermal 
density profile is always maximal in the center of the trap. 
In conclusion, we studied the effect of strong intcrac- 
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